Direct Measurements of Magnetic Twist in The Solar Corona 



A. Malanushenko 

Department of Physics, Montana State University 
<N: Bozeman, MT 59717, USA 

O 

(N : M. H. Yusuf 

(D; Berea College, Berea, KY 404O4, USA 

Ph 

^ . D.W. Longcope 

Department of Physics, Montana State University 
P^- Bozeman, MT 59717, USA 

ABSTRACT 

O 

^: 

^ ■ In the present work we study evolution of magnetic helicity in the solar corona. 

We compare the rate of change of a quantity related to the magnetic helicity in 
the corona to the flux of magnetic helicity through the photosphere and flnd 
that the two rates are similar. This gives observational evidence that helicity 
CN ■ flux across the photosphere is indeed what drives helicity changes in solar corona 

x/^ ■ during emergence. 

^ ■ For the purposes of estimating coronal helicity we neither assume a strictly 

^ ■ linear force-free fleld, nor attempt to construct a non-linear force-free fleld. For 

each coronal loop evident in Extreme Ultraviolet (EUV) we flnd a best-matching 
line of a linear force-free fleld and allow the twist parameter a to be different 
^ ■ for each line. This method was introduced and its applicability was discussed in 



Malanushenko et al.l fcoOQl ). 



The object of the study is emerging and rapidly rotating AR 9004 over about 
80 hours. As a proxy for coronal helicity we use the quantity averaged 
over many reconstructed lines of magnetic fleld. We argue that it is approximately 
proportional to "flux-normalized" helicity where H is helicity and $ is 

total enclosed magnetic flux of the active region. The time rate of change of such 
quantity in the corona is found to be about 0.021 rad/hr, which is compatible with 



the e stimates for the same region obtained using other methods (ILongcope et al. 



20071 ). who estimated the flux of normalized helicity of about 0.016 rad/hr. 
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Introduction 



Magnetic helicity is generally accepted to be an important quantity in understanding 
evolution of coronal magnetic fields and studying solar eruptions. It is approximately con- 
served when conductivity is high, as it is expected to be in solar corona. This sets up an 
important constraint on evolution of a magnetic field. It is believed that helicity in excess of 
a certain threshold could be responsible for triggering magneto hydrodynamics instabilities 
and thus be responsible for coronal mass ejections (CMEs) (see iDemoulinI 120071 . for further 
discussion) . 

Helicity is a function of magnetic field H = i^(B), defined as a volume integral of A • B 
(where A is the vector potential and B is the magnetic field), provided no lines of magnetic 
field leave this volume. For situations when field lines leave the volume (e.g., B • fi 7^ 
on the boundary), such as the solar corona, helicity is defined relative to some reference 
field that has the same normal component o f magnetic field at the boundary of the volume : 
Hreiative = H{B,Bref) = H (B) - H{Bref) dserger fc FieldlllQsi [pinn fc AntonsenI llQSsh . 
A potential magnetic field (V x B^ef = 0) has the minimum possible energy, so using it as 
a reference means a non-zero helicity demands some free energy. Relative helicity, defined 
this way, is approximately conserved under motions of plasma internal to the volume of 
integration. 

Direct measurements of helicity in the corona remain an extremely challenging problem. 
They are usually performed by extrapolating magnetic field into corona using photospheric 
magnetic field as a boundary condition and then estimating the helicity of this field. The 
extrapolations are often performed in such a way that lines of resulting magnetic field resem- 
ble observed coronal loops evident in extreme ultraviolet (EUV) or soft X-Rays (SXR). The 
popular choices of magnetic fields are li near force-free (or constant-a) field s (see Section [2] 



for description) confined to a box (e.g. iGreen et al 



2002 



Lim et al.l 120071 ) and non-linear 



force- free field (NLFFF) extrapolations (IRegnier et al.ll2005l ). Both methods remain imper- 
fect. The main drawbacks of the linear force- free field are that its is c learly wrong for active 
regions with field lines of clearly different twists ( IBurnette et al.ll2004j ). and it places restric- 
tions on a depending on the domain size. The second method has problems dealing with 
solving non-linear equations and wi th the use of ambigu ity-resolved vector magnetograms 
in a non force- free photosphere (see iDemoulin et al.lll997l . for further discussion). Applying 
diflFerent e xtrapolation method s to the same data was found to produce significantly diflFerent 
solutions ( IDeRosa et al.l 120091 ). 



A quantity that is easier to measure than helicity is the fiux of helicity through the 
photosphere. It can be shown that for changing relative helicity, changes of magnetic fiux 
at the boundary are far more effective than internal electric currents in the presence of high 
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connectivity (lBergerlll984[ ). This allows one to express the change in the magnetic helicity in 
th e corona thro ugh apparent motions of photospheric magnetic fe atures, as first suggested 
bv lChael koOlh and later developed by iDemoulin fc Bergerl (120031 ) . 



The theoretical prediction that the coronal helicity is injected from underneath the 
photosphere is not supported by much observational evidence. Only a few works compare 



integrated helicity flux through the photosphere to coronal helicity. IPevtsov et al.l ( 120031 ) 



have studied the evolution of coronal a and found that its evolution is consistent with the- 



oretical estimates of such for an emerging twisted flux tube. iBurnette et al.l (12004 ) have 
found a correlation between coronal a of a linear force-free field (chosen to visually match 
most of SXR loops) and a measured in the photo sphere using vec tor magnetograms (aver- 
aged, in some sense, over the whole active region). i Lim et al. ( 2007 ) have found that helicity 
injection through photosphere (obtained using local correlation tracking) is consistent with 
the change observed in the corona after taking account, approximately, of helicity carried 
away by CMEs; they assumes linear force- free field and a typical value for a helicity in CME. 
Comparisons of the change in coronal helicity with helicity of interplanetary magn etic clouds 
have demonstrated bot h a clear correspondence of the two (IMandrini et al.l 12004 ) and a lack 
of it (iGreen et al.ll2002l ) . iGeorgoulis et al.l (120091 ) have compared integrated helicity fiux with 
total helicity carried by CMEs (by multiplying a typical helicity for a CME by the amount 
of CMEs per studied period) and found an agreement between the two values. 

In the current work we study an emerging active region over a long period of time with 
high temporal cadence. We compute a rate with which a proxy for coronal helicity changes 
and find that this rate is consistent w ith rate of helici t y inje ction through the photosphere, 
reported for the same active region bv lLongcope et al.l (120071 ). This is the second time (after 
Park et al.ll2010l ) that the helicity change rate in the corona over a long period of time has 
been found to be consistent with the helicity fiux in the photosphere. 

We approximat e the s tate of the coronal field using a method recently proposed by 



Malanushenko et al 



(120091 ). It approximates coronal loops with lines of linear force- free 
field like the above mentioned works, however, it allows a to vary from line to line. That is, 
each of the coronal loops is approximated by a field line of a different linear force- free field. 
Strictly speaking, this approach is in correct, as a superpos i tion o f linear force-free fields will 
not resemble a force-free field at all. iMalanushenko et al.l (120091 ) argued that it might work 
for some cases relevant to the solar corona and have sup ported the reasoni ng with tests on 
analytical non-linear (or non-constant-a) force- free fields ( Low fc Loij 1990) with amount o f 



twist comparable to one typical to solar active regions as reported by lBurnette et al.l (120041 ) 



This method of determining properties of coronal loops does avoid some of the prob- 
lems faced by linear or non-linear extrapolations. We will henceforth refer to it as (Non)- 
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Linear Force- Free Field or shortly (N)LFFF. It introduces fewer uncertainties than a linear 
extrapolation, as a is allowed to vary in space. Also it is immune to some of the draw- 
backs plaguing non- linear extrapolations; in particular, it does not use vector magnetograms 
(w hich are obtained in the photosphere, where the force-free approximation is questionable. 



see 



Demoulin et al.l ( 119971 ) for further discussion) and does not attempt to solve non-linear 



equations. 

The paper is organized in the following manner. In Section [2] we review the methodology, 
describe the measured quantity and establish its relation to helicity. In Section [3] we describe 
the data to be used. In Section H] we summarize results and findings from the data analysis. 
Section [5] contains the discussion of the results achieved. The Appendix explains tangent 
plane projection and why it was chosen for the present work. 



The Method 



Most active region coronal magnetic fields are believed to be in a force-free state, 

V X B = aB, 



where a is a scalar of proportionality (e.g.. iNakagawa et al.lll97ll ). Extrapolations of coronal 
magnetic field usually imply solving this equation using photospheric data as a boundary 
condition. When Va = 0, Equation ([T]) reduces to a system of linear equations and thus is 
called linear or constant-a force-free field. When a varies in space. Equation ([T]) represents 
a non-linear system and the solution is called a non-linear of non-constant-a force-free field. 

In the present wor k we use linear f orce- f ree fields confined to half space, computed using 
Green's function from IChiu fc HiltonI (|l977n. Su ch fields are less popular than linear force- 



free fields confined to a box ( j Alissandr akisi 1 1 98 ll ) because it takes much more computational 
time to build them. We use them, however, because they do not impose restrictions on 
a based on the size of the computational domain and because the field lines are allowed 
to leave the computational domain. We believe this is a better representation of coronal 
magnetic fields. 

The result of the procedure of the loop fitting is a set of field lines of constant-a fieldfl, 
each of them is a best fit to an individual coronal loop for every visible loop. This gives the 



^Note that even for non- linear force- free fields a = const along every field line, of Va • B = 0. This 
result is obtained by taking the dive rgence of both sides o f Eq. ([T]) and using V • B = and the identity 
V • (V X B) = 0, so that B • Va = (|Priest fc Forbeslboool ). 
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following set of parameters for every coronal loop that was successfully fit: its a, shape r(/) 
and the magnetic field strength along the loop B(/). 

Note that this method is principally different from calculating maps of = Jzj^z using 
vector magnetograms. (N)LFFF uses coronal morphology to derive currents in the corona; 
these could, in principle, be traced down to the upper chromosphere level and provided 
enough coronal loops are taken into account, a map of could be then derived. In the 
previous work we ha ve shown that such maps yield acceptable results on synthetic data. 
Burnette et al.l ( 120041 ) have done similar analysis, but they estimated one a that would best 
fit all of the observed loop of a single active region (that is, assuming the field of an active 
region is linear). They compared the results obtained from vector magnetograms (a^) and 
the results of the single-a coronal fit a^est and reported the correlation between the two. 

In the current work we do not compute a helicity, but rather a quantity related to it, 
aiLi/2 (where Li is the length of i-th field line and ai is its a). Most classical definitions of 
helicity involve a volume integral, which requires knowledge of magnetic field on a grid or 
analytically in a volume. (N)LFFF does not provide such gridded data or even volume-filling 
data. However, aL/2 is closely related to helicity. For example, self-helicity of a uniformly 
twisted torus is i7 = Tw^^^ where Tw is a number of turns that a field line mak es per unit 
length and $ is a net magnetic fiux ( IBerger fc Field! 1 19841 : iMoffatt fc Riccalll992l ). and for a 
thin cylindrical uniformly- twisted fiux tub^ it can be shown that 2t{Tw = aL/2 over axial 
distance L (Aschwanden, 2006). 

It is not immediately clear that the same expression could be used for a more complex 
magnetic fiux configuration when thin fiux tube approximation is not applicable. However, 
Longcope fc Malanushenkd (120081 ) have demonstrated how additive self helicity (helicity of 
a field relative to a potential field confined to the same domain) is consistent with an em- 
pirical function a(L^)$^/47r with (L^) being the average length of a field line in the do- 
maiilfl and $ being the total magnetic fiux in the domain. This was shown for a case when 
the thin fiux tube approximation was clearly inapplicable: for a linear force-free field of a 
quadrupolar field confined t o a box, whose domain is the field connecting two polarities. Later 
Malanushenko et al.l ( 120091 ) suggested that "fiux-normalized" additive self helicity could be 



treated as a generalized twist, Ha/^^ = Twgen for an arbitrary magnetic configuration. That 



^Thin flux tube approximation usually refers to a structure with well-defined axis, with diameter small 
compared to its length, and with radius of curvature large compared to its diameter. 

^Here "domain" is a volume occupied by field lines connecting two given footpoints, so there is no magnetic 
fiux across the boundaries of such a volume, except at the footpoints. 
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leads to a conclusion that 



2TvTWgen = a{Li)/2, 



(2) 



i n a si milar manner as in thin flux tubes, at least in linear force- free flelds. iMalanushenko et al 
(I2OO9I ) demonstrated that Twgen is equal to Tw for a thin flux tube and that it behaves like 



Tw in arbitrary magnetic conflgurations, for example, serving as a kink instability threshold. 
We hereafter refer to Twgen = Ha/^^ as "twist". 



To the best of our knowledge, such relationship between a and helicity has not been 
established for general non-linear force- free flelds. However in a thin non-uniformly twisted 
cylinder it might be expected that Twgen would be proportional to axial length with a 
constant of proportionality, that should reduce to a if the cylinder were uniformly twisted. 
For a more complex conflguration Twgen might be expected to be proportional to some length 
scale times a constant that has dimensions of a. We choose to use a quantity 
averaged over many fleld lines as a proxy of Twgen and flnd that it changes consistently with 
the injection of coronal twist. 

We compare this qu a ntity to the normalized helicity flux across the photosphere , mea- 



sured by lLongcope et al.l (120071 ). They follow the results of I We 



sch fc Longcopel (120031 ) . who 



start with the commonly accepted expression for helicity flux ( 

dH 



Berger fc Field Il984l ) 



= 2 / [(Ap . B)v - (Ap . v)B] . MS, 



(3) 



here n is a normal vector pointing into the volume and Ap is assumed to be in Coulomb 
gauge with Ap • fi = 0. They split it into two terms, that depend on the normal and 
tangential components of the velocity on the boundary. These two terms correspond to the 
change in helicity due to the emergence of the new magnetic flux and to the surface motions 
of the existing magnetic flux respectively. Assuming the boundary is z = plane, the second 
term can then be written as 



dH 



2 / (Ap ■ Y^)B,dS. 



(4) 



Welsch fc Longcopd ( 120031 ) further concentrate on this term only, thus discarding the helicity 
flux due to change in the magnetic flux. They split the flux of unconflned self helicity (of 
a fleld conflned to half space relative to the potential fleld in half space) into two parts: 
the "spinning" and the "braiding" contributions, as illustrated on Figure [H The flrst one 
comes from rot a tion of footpoints about their axis and its change rate was expressed by 



Longcope et al.l ( 120071 ) through the average spinning rate 9 spin as: 



Tja 

spin 



2^ a spin") 



(5) 
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for a-th footpoint. The second one comes from the relative rotation of footpoints and its 
change rate could be expressed through the average tilt angle 9braiding- 



1 



nab 



braiding ^ ^ '^braiding' 



(6) 



hrn.id.inn 



are n ot equivalent to fluxes of "twist" and "writhe" helicities ( IBerger fc Field 



19841 : iMoflFatt fc Riccalll992l ). In addition to the discussions in the previously mentioned pa- 
pers, we would like to provide an illustrative example of this statement. Consider a thin 
untwisted half torus from Figured] (top row), deformed in the following way: flrst, each of 
the footpoints is rotated about its own center by 9 (second and third rows), then the foot- 
points are rotated about each other by 9 (bottom row). This transformation is equivalent to 
a rotation of the whole torus as a rigid body and adds no inner twist to it. However, an ob- 
server, who believes that Hi^^ aiding is equivalent to the flux of writhe helicity, might consider 
that the writhe helicity has increased by ^^^/tt, while the final configurati o n is t he same 
as the starting configuration in a rotated reference frame. iLongcope et al.l (120071 ) argued 
that "spinning" and "braiding" fiuxes of helicity might be produced by different physical 
processes. 



Combining Equations ([5j) and ([6j) and noticing that for a magnetic fiux balanced dipole 
= we get: 

(7) 



2^ spin ' ""spin "braiding ^ 

From Equations (Ej) and ([71) it follows, that if the quantity {aiLi) /2 is roughly proportional 
to Ha/^^ in a non-linear force-free field. 



its change rate should be proportional to 



27r 



(ttiLi)/2, 



(8) 



na _|_ nb 

^ spin ^ spin 



OfjCib 

braiding 



Note that Equation ([71) only represents the term given by Equation (|lj), which assumes 
constant magnetic fiux. Since we are interested in comparing results of the photospheric fiux 
of helicity with ours, we follow the same reas oning and neglect the effects due to the change 
in magnetic fiux as well. It was pointed out by lDemoulin fc Bergerl (120031 ) that the horizontal 
motion of fiux tubes would, in principle, include both contributions to Equation ([71), including 
that from vertical fiow (emergence). It therefore follows that if tracking methods used 
to measure horizontal speeds tracked magnetic footpoints, then their use in Equation ([71) 
would also capture bot h contributions. Te sts of these tracking methods, however, casts some 
doubt on this premise (jWelsch et al.ll2007l ). Local Correlation Tracking seems to capture the 
horizontal fiow speed far more accurately than it captures the contribution of vertical fiow. 
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Its application to Equation ([7j) thus captures most reliably the first term, neglecting, or 
underestimating, the contribution of emergence. 

{aiLi)/2 is really meant to represent Twgen and thus be relevant only to the additive 
self-helicity i/^, that is to helicity of the field in a subdomain relative to the potential 
field in the same subdomain, and Htotai is the helicity of the field in half space relative to 
potential field in half spa ce. However, the diflFerence H+ntni — HA only depends on the shape o f 
the subdomain and both iLongcope fc Malanushenkd ( 120081 ) and iMalanushenko et al.l ( 120091 ) 
have studied the cases with subdomains of relevant shapes and found this diflFerence to be 
small (except at the pre-eruption state). 



3. The Data 

The present study is devoted to AR 9004. It was emerging, as its magnetic flux was 
increasing. At the earliest time when the magnetograms were available it had about 30% 
of its final magnetic fiux. The other feature of this active region is that its footpoints were 
rotating about each other and about themselves at measurable pace. 

While we have only analyzed AR 9004, we included AR 9002 in the study. There are 
several reasons for this, such as proximity of the two active regions and the appearance of 
coronal loops interconnecting them. Also, while many loops could definitely be attributed 
only to one active region, to the other one, or connecting the two, the connectivity of 
other loops was not so clear. We included both ARs in the computational domain and 
reconstructed all coronal loops in TRACE 171 A field of view (that had both active regions). 
For the further analysis we only considered the loops that were found to start and end at 
AR 9004. 

We chose 21 full-disk MDI magnetograms in the time range between 2000-05-18 12:48 
and 2000-05-22 03:12. The first and the last images are shown on Figured 

The magnitude of the magnetic field on each magnetogram was corrected for the line-of- 
sight factor, assuming the magnetic field to be purely radial: ^^(x, y) = Bios{x^ y) / ^1 — p(x, yY / 
where p(x, y) = y^x^~+^ is the plane-of-the-sky distance from the disk center to the point 
(x, v). We did not account for possible inclination of the niagneti c field in active regions 
( HowardI 1991 ). The largest average inclination angle Howard ( 1991 ) found for growing ARs 



to be 25° for the leading polarity. The error in the strength of the photospheric magnetic 
field that we make assuming the field is purely radial is thus about cos 25° ^ 0.9. It is not 
clear how big of a diflFerence would it make, but in Malanushenko et. al. (2009) we report 
uncertainties of other origin that are bigger than 10%. We acknowledge that the analysis 
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in the present work has multiple sources of uncertainties, mentioned in this paragraph and 
further through the text; they could have a cumulative effect resulting in the real uncertainty 
bigger than 10%. 

In the next step, each magnet ogram was then remapped to the disk center in ortho- 
graphic projection (see Appendix). This was done to correct for the foreshortening distor- 
tions of the magnetogram. This step was necessary as during the observed time sequence the 
center of the couple 9002/9004 has traveled approximately from 15°N30°E to 15°N16°W. 

For each magnetogram we selected several (typically two) TRACE 171 A images taken 
within 30 minutes of the magnetogram. The selection criteria was the visibility of many 
distinct coronal loops. On each TRACE image as many loops as could be discerned were 
^^traced" , or visually approximated with a smooth curve, a two-segment Bezier spline (e.g.. 



Prautzsch et al.ll2002h . 



Each magnetogram on the tangent plane was rebinned on a coarser grid, half the 
size in each dimension, and than used to construct 41 linear force-free fields with a G 
[—0.05, 0.05] arcsec"-^ and with step of Aa = 0.0025 arcsec"-^ in a box with z^ax = 200 arcsec 
and Az ?^ 1.5 arcsec. We used a threshold of lOOG, that is, we did not account for weaker 
magnetic fields. The reason for that was to save computational time. The objects of interest 
were coronal loops of the scale of an active region, and such threshold is not likely to intro- 
duce significant disturbance in them. Both ranges of a and z have proven to be sufficient to 
reconstruct most of the loops. 

For every "traced" coronal loop and for many linear force-free fields we have computed 
many field lines along the line of sight that cross the midpoint of the loop. Figure [3] illustrates 
such field lines for one of the linear force- free fields. All of these field lines were discarded 
except for the "best-fit", shown on Figure [3] (right) in yellow. The selection procedure was 
the following: each field line (that is defined by a of the linear force-free field it belongs to 
and the coordinate along the line of sight h) was projected onto the plane of the sky. Then the 
mean distance rf(/i, a) was computed between this projection and the traced loop. Then the 
minimum of rf(/i, a) on a grid of given value s a an d h was found following a semi-automatic 



algorithm described in iMalanushenko et al.l (120091 ) . 



We found that for many loops (^55% of all 301 loops) the general shape of the rf(/i, a) 
parameter space matches one of the types described in the fitting scheme. As for the "un- 
described" types, we chose not to ignore them, but try to identify the non-hyperbolic valley 
and choose a local minimum on it that resulted in visually better fit. 

Among those "underscribed" types of parameter spaces there was one that was found to 
occur frequently enough (^ 36% of all loops) to warrant new classification. The "anomalous" 
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(or "non- hyperbolic" ) valley in the parameter space of this type looked like two branches at 
each side of a = line joined together (see Figured]). According to the above mentioned 
algorithm such a loop had to be either ignored or the local minimum on this valley that 
corresponded to lowest height was to be chosen. We found that such a choice results in a fit 
visually much worse than chosing the global minimum of this valley. We thus proceded with 
selecting the global minimum on the non-hyperbolic valley. 

This type of the parameter space had not been found on tests in analytic data and 
there is no clear indication that it is trustworthy, except for the better visual correspondence 
to the coronal loops. It might represent a distortion of "conventional" shapes due to the 
inclination of the line of sight. We would like to mention that such classification of param- 
eter spaces is subjective and sometimes ambiguous, but it has found to yield statistically 
reasonable results when tested on analytical fields. As for the new type of solutions, we also 
would like to mention that many previous studies relied on visual comparisons. For example. 



Burnette et al.l ( 120041 ) selected a linear force-free field with a^est that seemed to match visu- 
ally many coronal loops; they repeated same procedure for many active regions and found 
that a^est is generally consistent with the average a inferreded from vector magnet ograms. 



Green et al. 


(2002 


); 


Lim et al. 


(2007) 



as 

much attention. So while about ^55% of all our datapoints are backed up with tests on 
analytical fields, ^36% of them are obtained using methodology that has been previously 
used in many similar studies. 

This procedure resulted in a "best-fit" line of a linear force-free field for every coronal 
loop. The quality of the fit was visually judged and assigned a subjective grade of C (poor 
fit), B (good) or A (perfect match). For further analysis, only loops that had a fit quality B 
or A were used. For example, the fit on Figured] was given a quality grade B and the fit on 
Figure [3] was judged to be of quality A. 



4. Results 

We have performed the reconstruction procedure for the total of 303 coronal loops of 
AR 9004. 61 of them (20%) had quality B and 219 (73%) had quality A, so in total there 
were 280 (93%) successful reconstructions. The real percentage of the successful fits might 
be slightly lower, as some coronal loops seeming to belong to AR 9004 on EUV image were 
reconstructed as open field lines or as field lines interconnecting ARs 9002 and 9004, and 
some of those fits have failed as well. We could have estimated how many loops seemed to 
have incorrect reconstructed connectivity, but we chose not to do so. The reason is that it 
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would involve visual estimation of the connectivity that may or may not be right on its own. 
The following evolution was observed: 

• Up to t = 35.5 hrs. All times are in hours since 2000-05-18 00:00:00. Images 1-7 
on Fig. [8] the loops had mainly negative (left-handed) twist and the median twist has 
decreased in magnitude. 

• From t = 48.2 hr to t = 51.7 hr (images 10-11 on Fig. [8]) had twist of both signs^ 
with negative twist predominantly in the eastern half of the dipole and positive twist 
predominantly in the western half. The median was close to zero. (Only half of the 
region was seen on images 8 and 9, so we could not draw any conclusion about spatial 
structure of the twist from t = 38.0 hr to t = 41.3 hr.) 

• Most of the loops on image 12, t = 54.6 hr, had twist of positive sign. 



The lo ops on image 13, t = 58.9 hr, were poorly fit. As suggested in lMalanushenko et al. 



( 120091 ) , this might indicate that the field was strongly twisted or maybe strongly non- 



linear. 

Almost no loops of AR 9004 were observed and well fit on images 14-16 {t = 64.0 hr 
to t = 70.3 hr). It is possible that the existing loops were outside of TRACE field of 
view. 

Most of loops on images 17-21 {t = 72.2 hr to t = 99.4 hr) had twist of positive sign 
(or right-handed twist). The median twist appeared to be slightly increasing. 

The emergence of AR 9004 was evident from the evolution of the magnetic fiux. The 



magnetic fiux, according to lLongcope et al.l ( 120071 ). was steadily increasing until about t = 55 



hrs, and after that it exhibited a slight decrease. This decrease was concurrent with the drop 
in the combined twist angle on Fig. [5l 

For any given time there was a wide distribution of measured twist. There was, however, 
a general trend for the twist to increase. While most of the loops were negatively twisted at 
early times, most of loops at later time exhibit positive twist. 

We summarize the findings in Figure [61 Each point there represents aL/2 of a single 
loop versus time when it was measured. At every time we find a mean twist and mean of 
the absolute deviation (this way the time frames with many data points would not be given 
more weight in the fitting procedure than the time frames with few data points). We then fit 
a line to the twist as a function of time using least absolute deviation fit to the means using 
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the mean of the absolute deviation. We only fit the line in the time interval from t = 10 hrs 
to t = 55 hrs. The start time is the earliest time available and the final time is where the 
photospheric twist starts to decrease. We have measured aL/2 to increase at about 0.021 
rad/hr. 

We can obtain a sense of uncertainty in this measurement in several ways. Lacking a 
knowledge of either the magni tude or distributi on of the measurement errors we can employ 



the method of bootstrapping (iPress et al.lll986 ) which does not require such knowledge. We 



generate synthetic data sets by randomly resampling the actual data, permitting duplica- 
tion, and then repeatedly fit the synthetic data using LAD on the co-temporal means. The 
result is a distribution of measured rates whose central two quartiles are contained within 
(0.020, 0.035) rad/hr, as shown on Figure [3 A second estimate comes from the maximum 
likelihood formalism. A LAD fit is equivalent to the fit of maximum likelihood under the as- 
sumption of exponentially-distributed, additive errors with identical mean deviation equal to 
that of the best fit: 0.237 rad. The central 50% of the likelihood distribution is encompassed 
by an (asymmetric) range of slopes (0.017, 0.035) rad/hr, similar to that from bootstrapping. 
Finally, we attempted several different methods of fitting linear trends to the data including 
least-squares fit, fits to all data points, to the medians of each time (rather than the means), 
and using only A-quality points. Each method yielded a diflFerent slope, but the collection 
fell within the range (0.016, 0.038) rad/hr. 



Figure [5] shows the twist injection rate as measured by iLongcope et al.l (120071 ) as a 
black curve. It also shows individual spinning (red for positive and blue for negative) and 
tilt (green) angles. We fit a line to the combined angle using least absolute deviations to 
the same time interval as on the coronal measurements and find a twist injection of about 
0.0160 rad/hr (what corresponds to the -0.0160 rad/hr change in 9^^-^ + 9~p-^ — "^.d^raiding^ 
Equation [Tj) . 

We conclude, that the rate of change of the coronal twist is consistent with the twist 
injection through the photosphere. 



5. Discussion 

In the current work we have observed how helicity fiows into the corona through photo- 
spheric motions. In the beginning of the time sequence coronal loops of AR 9004 appear to 
have negative helicity and after about 60 hours all the coronal loops appear to have positive 
helicity. We have observed that when a negatively twisted field is subject to the injection of 
helicity of positive sign, the average twist increases and eventually passes through zero (like 
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at t G [40,50] hrs). But this does not mean the magnetic field becomes potential. Rather, 
it passes through a complex non-linear equilibrium that has parts with distinct positive and 
parts with distinct negative twist. 

There were two weak fiares associated with AR 9004 in the examined time interval, 



according to Solar Monitor (IGallagher et al.ll2002r). Th ere was one C2.7 fiare at t = 52.6 hrs 



and then a CI. 9 at t = 88.6 hrs. iKusano et al.l (120031 ) have proposed a model in which the 
field of an active region, featuring twist of both signs, undergoes magnetic reconnection. This 
process is accompanied by a fiare and results in an untwisted magnetic field (that is, with 
zero helicity). The fiares in AR 9004 do not seem to be associated with reliably measured 
drops in twist. Moreover, the AR that we study does have twist of both signs at t G [20, 55] 
hrs and yet does not seem to relax to a potential state at later times, after t = 72 hrs. 

It is worth noting, however, that in the period of t G [58.9, 70, 3] hrs there were almost 
no loops that were successfully fit. This opens the question of what might have happened 
at that time. It could be argued that the field indeed had relaxed to the potential state and 
its further positive (right-handed) twist was injected through the photosphere, but at least 
the apparent injection of negative helicity within about t G [55, 70] hrs (evident on Figure [5j) 
suggests that this might not be the case. 

We have studied the time rate of change of the following quantity: (a^L^/2) averaged 
over many reconstructed field lines. In Section [2] we argued that it might be proportional 
to the additive self helicity of a non-linear force- free field in the similar manner as aL/2 is 
proportional to twist helicity of a thin fiux tube and as a{Li)/2 is related to additive self 
helicity for a linear force- free field. 

We have found that the time rate of change of (a^L^/2) (and arguably a generalized 
twist, Twgen = Ha/^^) is found to be about 0.021 rad/hr. This rate is similar to the time 
rate of change of the fiux-normalized total helicity Htotai/^^ that was found to be 0.016 
rad/hr. 

The diflFerence between Ha and Htotai might be responsible for the fact that even though 
Htotai starts to decrease after about t = 55hrs, the twist derived from coronal loops remains 
of a positive sign after t = 70 hrs (or arguably t = 55 hrs, as explained below) and does not 
seem to obviously decrease. We have noticed that a lot of bright coronal loops connecting 
ARs 9004 and 9002 appear at about t = 75hrs and later at about t = lOOhrs. It is possible 
that while photospheric helicity injection changes sign, magnetic reconnection that happens 
afterward changes the balance between "twist" and "writhe" helicities, in the sense of Ha 
as a "generalized twist helicity" and Htotai — Ha as a "generalized writhe helicity". In 
the current paradigm magnetic reconnection results in a decrease of magnetic energy. It is 
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also true that linear force-free fields with lar ger a are typically thought of as having larger 



magnetic energy than those with smaller a ( jAlylll992l ). However the field of AR 9004 was 



found to be non-linear, which means , for the same system there exists a linear force- free 



field with lower energy ( lWoltjerlll958l ). A possible scenario for the observed phenomena is 
the following: magnetic reconnection between ARs 9004 and 9002 lowers the total magnetic 
energy and possibly turns AR 9004 to a linear (or nearly linear) force-free field after about 
t = 70 hrs. Within the time interval approximately t G [55, 70] hrs the total helicity of AR 
9004 decreases due to photospheric motions (see Figure [5j), however, the general trend of the 
self-helicity within t G [55, 100] hrs does not seem to demonstrate a similar decrease, except 
for maybe t G [60, 70] hrs, but the data in that range are poorly sampled; each greyscale bin 
has only a few data points and the means are computed for individual time frames that have 
even fewer points. The means seems to decrease, however, the spread in values is similar to 
that in t G [50, 60] hrs. Our own belief is that these two points are less reliable than the rest 
of the data. This is consistent with a possible interpretation of the data provided earlier in 
the text. This decrease in the total helicity possibly signifies a decrease a "generalized writhe 
helicity" component. The latter appears to be only a function of the shape of the domain 
containing field of AR 9004. So both magnetic reconnection (that has decreased the volume 
of this domain by reconnecting some of the magnetic fiux to AR 9002) and photospheric 
motions might have contributed to that. 

This work was supported by NASA under grant NNX07AI01G and NSF under award 
ATM-0552958. We are grateful to Carolus Schrijver for the name "(non)-linear force-free 
fit" and the anonymous referee for the useful suggestion that led do a significant improve of 
the paper. 
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Fig. 1. — A schematic illustration of how braiding and spinning motions of footpoints contribute 
to helicity. (top row) An untwisted half torus, (two rows in the middle) Each of the footpoints is 
rotated about its own center by angle — 9- — 9 (spinning motion), (bottom row) If then the 
footpoints are rotated about each other (braiding) by 9i)j. — 9^ the resulting configuration would 
have no helicity in it, in agreement with Equation [71 
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Fig. 2.— We studied the region with AR's 9002 and 9004 from 2000-05-18 12:48:03 till 2000-05-22 
03:12:03. AR 9002 was an old diffuse active region and AR 9004 was emerging and rapidly rotating 
active region. 
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Fig. 3. — (left) - TRACE image in the plane of the sky and de-rotated MDI images shown with 
their respective angles. For every coronal loop (yellow) the algorithm browses through several 
different linear force- free fields. In each field it traces field lines (cyan) along the line or sight 
(red). The field lines are than projected back to TRACE image and compared to the loop using 
the semi-automatic selecting algorithm. Best-fit field line is shown in darker shade of blue. (Field 
lines are of a = 0.025 arcsec"^, that corresponds to the best-fit a for this particular loop, so the 
dark blue line is the resulting best-fit line for this loop.) Gray box shows the actual computational 
domain, (right) - a fragment of the TRACE image in the plane of the sky, the traced loop (yellow), 
the projection of the best fit to the plane of the sky (cyan) and the midpoint of the loop (red). 
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Fig. 4. — A type of parameter space, not described in iMalanushenko et al.l (l2009l ) on tests with 
analytic fields, (a) — The parameter space. White diamonds: local minima in columns. White 
diamonds with red dots: selected "non-hyperbolic" valley. Three colored squares: points that 
correspond to the local minima on this valley. Blue dotted lines: hyperbolas h — mi / n — 
0, =bl, ±2, .... (b) — d{a^h) along the selected valley, (c) — A fragment of TRACE 17lA image 
with the traced loop (dashed red) and three field lines, corresponding to the three local minima. 
It seems that the one, corresponding to the lowest h (magenta) is visually a much worse fit than 
the one, corresponding to the global minimum on this valley (cyan). We consider the valley to the 
left of the 'hump' a different valley, that could have been a degenerate hyperbolic. Based on this, 
we decide that the cyan line should be used. We update the fitting algorithm with the parameter 
space of this type. 




Fig. 5. — Twist rate, as measured by lLongcope et alJ (120071 ). The data points within time interval 
from to = 12 hours to ti = 57 hours were fit to a line with least absolute deviations. The beginning 
of the time interval was chosen to be the lower boundary for (n)lfff data and the end of the time 
interval was chosen where the linear twist injection was no longer obvious. 
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Fig. 6. — Twist of coronal loops versus time, using (n)lfff reconstruction from 
(120091 ). (left) — Diamonds show twist of individual loops (larger correspond to quality A and 
smaller to quality B). Blue squares are medians for each individual time. The line shows least 
absolute deviation to the means with means of the absolute deviations fit to the diamonds within 
the selected time range (same range as on Figure Ej). (right) — A histogram showing time-twist 
distribution of coronal loops makes the trend evident. Black color corresponds to 20 points or more 
and white corresponds to one or no points. 
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Fig. 7. — Evaluation of the uncertainty of the twist injection rate using bootstrapping method. 
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Appendix: Tangent Plane Projection 

A tangent plane projection is an orthographic projection onto a plane, tangent to the 
Sun at the point which is called the center of the projection and the projection of the 
solar North directly upwards. It could easily be obtained from the plane-of-sky using the 
following transformations. I n the starting pla ne-of-sky the coordinate system is assumed 



to be heliographic-cartesian ( IThompsonl l2006l ) . That is, Cartesian, with the origin at the 



Sun's center, x axis directed towards solar West in the plane of sky, y axis directed towards 
solar North in the plane of sky and z axis directed towards the observer. The image plane 
is then rotated as described by Equation [91 first, by —60 (where 60 is solar B-angle) about 
X-axis, then by —0c (where (/)c is the longitude of the desired projection center) about the 
new y-axis, then by —9c (where 9c is the latitude of the desired projection center) about 
the new x-axis. After this sequence of rotations solar North would lie on the new ^-axis 
and solar West would lie on the new x-axis. The last step is to perform the orthographic 
projection in the new z direction, that is, simply to set z coordinate to for all points in 
the visible hemisphere. 



AO \ 

A^{e) = cos(^) -sin(^) 

VO sin(^) cos(^) / 

(cos(0) sin(^)\ 

10 

— sin((/)) cos((/))/ 



A = A,ibo)A,i-0Myi-4>c)AA-bo) (9) 

This projection is non-conformal (it does distort local shape) and distorts the distances. 
But those distortions are independent of the location of the point of tangency and only 
increase with the size of the desired box. For example (see Figure [9j), if one considers a 
point on the sphere, which radius-vector from Sun's center r makes an angle 7 with the 
radius- vector from Sun's center to the projection center, Fq, then the distance on the sphere 
between r and Fq is Rj and the distance between the projections of these two points on the 
tangent plane is i?sin7. For 7 = 30° (or the box of 60°, which is enough to fit most active 
regions), the foreshortening factor would be about 5%. 

For purposes of tracking long time sequences and dealing with off-center solar images, 
heliographic coordinates (with longitude (j) along horizontal axis and latitude 9 along vertical) 
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Fig. 9. — Distortions of distance on tangent plane projection are small if the size of the extracted 
region is small. For example, if a distance between the point of tangency and another point on the 
Sun is i?7, then (if 7 < 7r/2) the distance between the center of the projection and the projection 
of that point is i?sin7. 



are o ften used, with the longitude either equally spaced along the verti cal axis (iThompsgn 
2006 ) or with spacing that changes with the distance from the equator (Welsch et al. 12009 ). 



The first one is commonly referred in cartography as plate Carree projection and a p articular 



case of the second one frequently used in Solar Physics is called Mercator projection ([Pearson 



1990l ). 



A tangent plane projection, such as we use, is somewhat less common. The first reason 
why we chose it is that the distortions are independent of latitude (unlike for Mercator 
or plate Carree projections, which have systematic latitude-dependent errors), as explained 
above. The second reason why we choose the tangent plane projection is that it is distorts 
global shapes less than cylindrical projections. Since the shape of a coronal loop is crucial 
for determining its a, we believe it is more important to preserve shapes globally rather 
than locally, and orthographic projections are better in this sense than cylindrical ones. For 
example, let us consider a 30° x 30° "square" ABCD^ centered at 0°W, 15°N (roughly the 
size and position of ARs 9002 and 9004 when they pass through the central meridian), as 
shown on Figure [lOl The lengths of the arches on the sphere are AB = EC = AD = 
fi?0 ^ O.524i?0 and CD = fcosfi?© ^ 0.454/?©, so CD/AB ^ 0.86. In a cylindrical 
projection CD/AB = 1 and on a tangent plane centered at 0°W, 15°N, CD/AB ^ 0.87 
(the latter is obtained by writing parametric equations for the projections of the box's sides 
and computing their length using standard methods). The foreshortening factors on the 
tangent plane are about the same for all four sides and are about 0.99. We thus believe that 
by using tangent plane projection in this particular case we make about 1% error due to 
the projection effect, as opposed to about 15% error in a cylindrical projection (given for 
CD/AB] the ratio AB/ BC would differ depending on the particular type of the cylindrical 
projection). 
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Fig. 10. — Illustration of heliocentric-cartesian (orthographic) and heliographic projections. 
ABCD is a 30° x 30° "square", centered at 0°W, 15°N (roughly the size and position of ARs 
9002 and 9004 when they pass through the central meridian) in three different projections. (Left) 
- in the plane of the sky (orthographic projection, centered at 0°W, 0°N and neglecting the &-angle 
for illustrative purposes), (middle) - in the tangent plane (orthographic projection, centered in the 
middle of ABCD, i.e., 0°W, 15°N), (right) - in heliographic coordinates (plate Carree projection). 
Dashed lines are lines of constant latitude and longitude. 
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